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Abstract 


We propose a new approach for solving nonlinear Klein-—Gordon and 
sine-Gordon equations based on radial basis function-pseudospectral method 
(RBF-PS). The proposed numerical method is based on quasi-interpolation 
of radial basis function differentiation matrices for the discretization of spa- 


tial derivatives combined with Runge-Kutta time stepping method in order 
to deal with the temporal part of the problem. The method does not require 
any linearization technique; in addition, a new technique is introduced to 


force approximations to satisfy exactly the boundary conditions. The intro- 
duced scheme is tested for a number of one- and two-dimensional nonlinear 
problems. Numerical results and comparisons with reported results in the lit- 


erature are given to validate the presented 


method, and the reported results 


show the applicability and versatility of the proposed method. 
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1 Introduction 


The Klein—Gordon equations play an important role in mathematical physics, 
for example, plasma physics, solid state physics, fluid dynamics, and chemical 
kinetics [12, 14, 4]. The sine-Gordon equations arise in many different appli- 
cations and model a soliton wave, which occurs in several physical situations 
[14]. The theory of plane or cylindrical Josephson junction requires the so- 
lutions of sine-Gordon equation in two dimension space [11]. The Josephson 
junction model, which has various applications in physics, electronics, and 
so on, consists of two superconducting layers with a tunnel junction between 
them. This model can be described by the two-dimensional sine-Gordon 
equation. Here, we consider the following time dependent nonlinear partial 
differential equation: 


ur +adut+N(u) = f(z,t), (#,t)€ xR, (1) 


with initial conditions u(z,0) = uo(x), OQu(x,0) = vo(x), and Dirichlet 
boundary conditions such that (uo,vo) belongs the energy space H?(Q x 
R*) x L?(Q x Rt), where O C R4,(d = 1,2). We consider the above Cauchy 
problem associated to the semilinear wave equation, where a < 0 is a known 
constant, A = ear Or is the standard space Laplace operator, and NV 
stands for a nonlinear function. In the well-known sine-Gordon equation, the 
nonlinear function is given by V’(u) = sin(w). On the other hand, the higher 
order Klein—Gordon equation with quadratic and cubic nonlinearity is given 
by N(u) = Bu+yu* for k = 2,3. Here, we consider the numerical approxima- 
tion of one- and two-dimensional sine-Gordon equation and the quadratic and 
cubic nonlinear Klein—Gordon equations. The well-posedness of the Cauchy 
problem associated with the semi-linear wave equation has been proved in 
(6, 29]. Many numerical methods such as the finite-difference method, finite 
element method, spectral method, and so on, have been employed to solve 
nonlinear Klein—Gordon equations mentioned above[16, 44, 47]. 


Recently, several numerical methods have been proposed based on radial 
basis functions (RBFs) in order to solve various types of ordinary and partial 
differential equations|2, 3, 1, 34, 28, 40, 39, 10, 27, 30, 23, 15, 36, 5]. Usually, 
RBF methods need to solve the ill-conditioned linear system in order to solve 
partial differential equations (PDEs). To overcome this difficulty, various 
methods have been proposed, such as the multiquadric quasi-interpolation 
(MQQI) method [32]. Here, we will study the numerical solutions of the 
nonlinear Klein—Gordon and the sine-Gordon equations based on a combina- 
tion of radial basis functions (RBF) and pseudospectral method MQQI. The 
differentiation matrices based on MQQI method, are constructed to discretize 
the spatial derivatives, and Runge-Kutta time stepping method is used to 
deal with the temporal part of the problem. The method does not need to 
use any linearization technique. In addition, a new technique is introduced 
to force approximations to satisfy exactly the boundary conditions. The pre- 
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sented method is tested for one- and two-dimensional nonlinear problems in 
order to show its efficiency and accuracy. Numerical results and comparisons 
with reported results in the literature are given to validate the presented 
method. 


2 Multiquadric Quasi-interpolation 


Quasi-interpolation methods are very useful tools in the analysis of scat- 
tered data both in theory and in application areas such as the computer 
science, economy, medicine, and geology. Multiquadric functions were intro- 
duced by Hardy [25], and then Franke [18] used them in numerical exper- 
iments. Authors in [35, 8, 7, 22, 19, 20, 21] successively proposed a num- 
ber of quasi-interpolation schemes and discussed their convergence. Useful 
quasi-interpolation operator Lp f was proposed by Wu and Schaback [46]. 
They discussed the convergence and shape preserving properties of this oper- 
ator. They proved a convergence theorem for interpolated functions f € C?. 
Zhang, Liang, and Li [48] proved the convergence of quasi-interpolation op- 
erator £pf for a wider range of functions. 


Assume that there are N points in the bounded interval 2 = |a,}] as 
follows: 
A=%,<%<-+++<Uay_-1<4¢y =). 


Let 
h=max{(x;-—2;-1), j =2,...,N}. 


In [46], the MQ quasi-interpolation operator £p has been defined: 


(Lof)(«x = Ye) (2) 
where 
n@)= b+ te a 
wo) esa 
$;(z) = ey ( scam ? oe serge iv = 2s (3) 
xz—xr)—y x wWn—1(©)—-wn_—2(x 
n—1(8) = “ary =ay—4) caren as 


bn (a) 1, n-1(x%)—-(t@n—-2) 


— 2 2(a@n—-2%N-1) 


and w,;(x) are the multiquadrics RBF (MQ-RBF) as follows: 


? 


W(x) =4/e+(a-—24,;)*, 7 =2,...,N—-1. (4) 
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The univariate MQ quasi-interpolation scheme can be extended into higher 
dimensions simply by applying it to each of the coordinates independently. 
In this paper, the scheme of dimension-splitting MQ (DS-MQ) proposed in 
[31, 33] for approximation of f(r), = (€,7) € R? is used, and this function 
will be approximated as follows: 


Ni, No N 
f(x) = f(E.n) = (Lof)(E.n) = S26 fisbi(©) G5 (n) = So fee (x), (6) 
i=1 j=l k=1 


where N = N, x No is the number of data points, fi; is the value of f at 
the intersection of the ith horizontal grid line and the jth vertical grid line, 
and ¢;(€) and ¢,(7) are defined by (3). Functions ®;,(2) = ®,(€,7) in (5) 
are obtained by product of ¢;(€) into $;(7) and f, = f(ae) = f(&, mK), & = 
1 N. 


pete 


3 Derivative approximation by the RBF-PS method 


The differentiation matrices have proved to be a very useful tool in the numer- 
ical solutions of differential equations. In this work, differentiation matrices 
are derived from the radial-basis pseudospectral (PS) method. PS methods 
are numerical methods with high precision, which use polynomials as ba- 
sis functions. Chebyshev pseudospectral method (CHEB-PS) is one type of 
them, which uses Chebyshev polynomials as basis functions. An important 
attribute of PS methods is that approximate solutions are obtained on a dis- 
crete set of grid points z;,7 = 1,...,N and in solving differential equations 
with PS methods, we will need to represent values of derivatives of u at grid 
points. The radial basis function-pseudospectral (RBF-PS) method involves 
the common approach of all weighted residual methods, that is, approxi- 
mating the field using a combination of some pre-defined functions termed 
as shape functions and some unknown coefficients, which are determined de- 
pending upon the governing equation of the problem. In the RBF-PS method, 
shape functions are computed by using RBF's. Fasshauer [17] proposed the 
application of infinitely smooth RBFs in pseudospectral formulation and in- 
terpreted the combined approach as a meshless RBF-PS method. 

For a given set of centers 7;,7 = 1,...,N, an MQQI approach of the 
solution u(x,t) of a PDE is given by 


€5 (1); (x), (6) 


un(a,t) = 


Me 


j=l 


where ¢;(x),j7 =1,...,N are defined by (3) for one-dimensional and ®,(z, y) 
are defined by (5) for two-dimensional. 
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Now, if we evaluate the approximate solution uy(z,t) at the grid points 
xj,i=1,...,N, and at fixed time ft, then we get 


N 

j=1 
Which can be written in the following matrix form: 

u(t) = A&(t), (8) 
where 
&(t) = [1 (t), tee En ()]’, 
Aij = 0;(2i), 459 = Vee cag dN, 

and 


u(t) = [un (21,t),... ,un(2n,t)]’. 


Suppose that £ is a linear differential operator relative to the x. By 
applying the operator £ to expansion (6) and use linearity of £, we have 


N 
Lun (x,t) = Db) E5(t)£9;(2). (9) 
j=l 
If we again evaluate (9) at the grid points 7;,i = 1,...,.N, and using matrix 
notation, then we get 
Lu= Ag&(t), (10) 
where €(t) is as above, Lu = [Luy (21, t),..., Lun (rn, t)]7, and the matrix 


Ag has entries £¢;(x;). Then, we use (8) to solve the coefficient vectors 
€(t) = A~tu, and then (10) yields 


Lu = AgA'u, (11) 


so that the operational matrix £ corresponding to linear operator £ is given 


by 
f=AgA*, (12) 


In order to obtain the differentiation matrix £, we need to ensure invertibility 
of the evaluation matrix A. This generally depends both on the basis func- 
tions chosen as well as the location of the grid points x;. For positive definite 
RBFs, it is well known that the evaluation matrix is invertible for any set of 
distinct grid points [45, 38]. For positive semidefinite radial basis functions, 
it is possible to construct a counterexample that leads to a singular matrix, 
though, such counterexamples are found rather seldom and with difficulty 
[45]. 
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Remark. In the special case if £ = A, then 

Au = De.u, (13) 
where 

De=AgA?, (14) 


x Ps - 2 
such that in one-dimensional A = ay and 


. a é 2 2 
Also, in two-dimensional, we have A = a =f a and 


PO; (x,y) , 0°; (2,y) 
Oxr2 I Oy? ) |(e,y)=(as yi)? 


(Ag)ij = ( 


4 New numerical scheme 


Now, we present the numerical scheme for solving the Klein—Gordon equation 
(17)-(20) by using the MQ quasi-interpolation operator Lp. In our scheme, 
the collocation scheme is applied at grid points x;, 7 = 1,...,N. Then, 
the highest order derivatives (second order in this paper) of the solution 
function are approximated by the differentiation matrix. Then, the fourth 
order Runge-Kutta method is employed to the approximate solution of the 
first order ordinary differential equations (ODEs) system. 


We first discuss our scheme for the nonlinear Klein—Gordon equation with 
quadratic and cubic nonlinearity in one-dimensional and two-dimensional 
space, and finally, we apply the numerical algorithm to solve the sine-Gordon 
equation. In the following, we consider the following nonlinear Klein—Gordon 
equation: 


Uee( 2, t) = aAu(a, t) = Bu(z, t) + y(u(a, t))F = FG; t); (17) 
for a < 0 with the initial conditions: 
u(z,0) =uo(z), «EQ, (18) 


uz(z,0) =vuo(z), 2EQ, (19) 


and Dirichlet boundary conditions: 


u(x,t) =g9(t), te(0,T], «xe dQ, (20) 
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where 2 is a bounded domain in R¢(d = 1,2) and 0 is the boundary of 
Q. We transform the problem (17)—(20) into a system of coupled first order 
equation in time given by 


u(x,t) = v(z,t), 


u(x,t) = —aAu(z,t) — Bu — yu* + f(a,t), ot 
with the initial conditions 
u(z,0)=uo(r), «rE, (22) 
u(x, 0) =v (x2), we, (23) 
and boundary condition 
u(z,t)= g(t), t€(0,T], «eon. (24) 


In the following, we give the structure of this scheme. 


4.1 Step one: Semidiscretization in space 


If we evaluate the MQQI approach of u(z,t) and v(z,t) at grid points 2, 
i=1,...,N and at fixed time ¢, then we get 


un(xi,t) = d- &(t); (ai), i=1,...,N, (25) 
j=l 
N 
un (zi, t) =>" Bi(t)o; (2), t=1,...,N. (26) 
j=l 
If we evaluate (21) at each node z;,7 = 1,..., N, and use (13) to approximate 


Au(a;,t) and equations display in a matrix form, then we get 


a = V(X, t), 
(27) 


W = —aDe.U(X,t) — BU(X, t) — yU*(X, t) + F(X, 0), 


where X = [1,%2,...,2n]?, U(X,t) = [ur(t),...,un()]?, V(X,t) = 
[vi(t),...,vw(t)]? and f(X,t) = [f(x1,t),..., f(zn,t)]*, in which u,(t) and 
v;(t) are approximations of u(z;,t) and v(2;,t),i = 1,...,N and Dg is the 
differentiation matrix defined in (14) and unchanged in time steps. In a more 
compact form, (27) can be written as 
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a = F,(t,U,V), 


ay = F(t, U, V), 


with initial values 
U(0) = u(X,0) = uo(X) 


V(0) = uz(x, 0) = vo(X), 
where 


F, ie U, V) = V(X, t), 


F,(t,U,V) = —aDg.U(X, t) — BU(X,t) — yU"(X, t) + F(X, 0). 


(28) 


(29) 


(30) 


Remark. In the above, the vector U* equals U QU ®--:@U and the 
symbol ®) denotes the component by component multiplication of two vec- 
tors. We can easily repeat the above step for nonlinear sine-Gordon equation 


by replacing the vector U* with the vector sin(U). 


4.2 Step two: Time discretization 


In this section, we apply fourth-order Runge-Kutta scheme to solve the sys- 


tem of ODEs (28). 


Let 
N 


um) =U = (un (entm) a = Oo eed) Ay 


j=l 


and let 
N 
VO) = VM = wy (ei tm) ha = O_ Bf os(ed) Ra 
j=l 


where eo =e stew) ae = 8; (tm), and tm = mdt. 
Then, from (28) and fourth-order Runge-Kutta method, we have 


pe) =U™M41(Ky, + 2K 91 + 2K 31+ K41), 


yimt)) _ yim) +4 1(Ki.2 + 2Ko2+2K32+ K42), 


where 

Kin = dt Fy (t™, ue. yr) , 
Koi _ dt F(t oe tt u™ tie . 
K3, = dtF y(t + 4,0 + $Ko1,V™ + $K22), 


(31) 
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Ka, = dtF (t™ + dt,U™ + K31,V™ + K30), 
and K;,9 is similar to K;,1,i1=1,...,4. 


4.3 Step three: Apply the boundary conditions 


The other half of the technique is the imposition of the boundary conditions. 
There are two basic ways of dealing with boundary conditions. The first is to 
use the basis functions that satisfy the boundary condition exactly. This is 
possible when the domain is regular and the boundary conditions are simple 
such as homogencous or periodic boundary conditions. In other cases, the 
second approach may be used, where some of the additional equations are 
used to enforce the boundary conditions. Here, we use the algorithm given in 
[17, 42] with some modifications. The boundary conditions will be imposed 
by modifying vector U("+) in (31). Replace those components of uit) 
corresponding to grid points at boundary domain with the exact boundary 
conditions at time level (m+ 1)dt. So we modified U“"*+”) in (31) with the 
following: 

ued) = Tg) 4 Ue) (32) 


where Ig is the matrix resulting from the replacement the rows corresponding 
to the boundary points with zero in the identity matrix and ugett) is the 
vector resulting from the replacement the components corresponding to the 
boundary points with the exact boundary conditions in the zero vector. 


5 Numerical examples 


In this section, several numerical examples are given to show the efficiency 
of our proposed method for approximating the solution of nonlinear Klein— 
Gordon equation and Sine-Gordon equation by comparing with other meth- 
ods. 

we first consider nonlinear Klein—Gordon equation (17)—(20), and finally, 
this scheme is applied for Sine—Gordon equation in both one-dimensional and 
two-dimensional cases. 

We use the univariate and two-dimensional MQQI Lp directly to approx- 
imate functions u™ and Au”. For all computations in this paper, we choose 
fixed values of the shape parameters of the kernel functions, and the rela- 
tionship between the shape parameter and the error norm is shown in the 
examples. In order to test the performance of the method and accuracy of 
approximate solutions, L.., £2, and root-mean-square (RMS) errors will be 
used, which are defined as 
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Loo error = ey |u(t5,tm) — un (%;,tm)|, 


N 
Lz error = h> (alate) —un(2j,tm))?, 
j=l 
Yijun (u(y, tm) — un (@4, tra)? 
RMS error = _ N ; 


for 0 <t», <7, where u and uy denote the exact and numerical solutions of 
the problem, respectively. 


5.1 Klein—Gordon equation 


In the case of one-dimensional nonlinear Klein—Gordon equations, we com- 
pare our numerical results with analytical solutions, solutions in [13], and 
those in [37]. Dehghan and Shokri [13] used thin plate splines (TPS), RBF 
collocation method. To approximate the solution function and second order 
derivative of the solution, Sarboland and Aminataei [37] used the integrated 
form of MQ quasi-interpolation scheme(IMQQI). They showed that their 
scheme requires fewer nodes to attain accuracy of TPS method, and in this 
work, it is shown that for presented examples, the errors in our method are 
better than those in [37]. In two-dimensional Klein—Gordon equations, we 
compare our numerical solution with the exact solutions, approximate solu- 
tion in [41], and we show that the numerical results are in good agreement 
with the exact solution and that our scheme requires fewer nodes to attain 
accuracy of MK-MLPG method. 


Example 1. Consider problem (17)—(20) with interval x € [-1,1], a = 
-1,2 =0,y=1, and k = 2. The initial conditions are given by 


ee eae 


and f(z,t) = —xcost + x? cost in the interval -1 < 2 < 1. The exact 
solution is given by 
u(x,t) = xcost. (33) 


The Dirichlet boundary conditions are extracted from the exact solution eas- 
ily. 

Table 1 presents L,, errors, D2 errors, and RMS errors in Example 1. The 
results are obtained by using N = 10 grid points, dt = 0.001, and c = 0.121. 
We compare our results with the reported results in [37]. As can be seen from 
the tables, the numerical results are in good agreement with the exact solu- 
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Table 1: Error norms of Example 1 at different times T using N = 10 grid points, 
dt = 0.001, and c= 0.121 
T Loo Lo RMS L~[37] L2[37] RMS[37] 
1 2.19755 x 10- = 2.1570 x 10-™ 1.44696 x 10-1 1.2590 x 10-5 2.0694 x 10-> 6.2397 x 10-% 
3 1.30590 x 107!4 1.23829 x 10-!4 8.30668 x 10-15 1.5428 x 10-5 3.7065 x 10-5 —-:1.1175 x 10~° 
5 2.72352 x 10-14 2.57676 x 10-14: 1.72854 x 10-4 3.3625 x 10-5 6.9684 x 107° 2.1010 x 10-° 
7 2.08028 x 107!4 2.03321 x 10-!4 1.36392 x 10-!4 3.7412 x 10-5 8.1943 x 1075 —- 2.4706 x 107° 


2.5x10-"4! 
5 2,x10-*4| 
a =1 
5 @ 15x10" == 
z —t=3 
[=] 14} 
8 LEW } — tz5 
5.x 10-5 —<ital 
0 
i -1.0 
* x 


Figure 1: Numerical and exact solutions of Example 1 (left) and error of method (right) 
with N = 10, dt = 0.001 


tion.The graphs of the numerical and analytical solutions for different times 
with dx = 0.2,dt = 0.001 in the interval —1 < x < 1 and using MQQI-PS 
(with c = 0.121) are shown in Figure 1. 


Example 2. In this example, we consider the one-dimensional nonlin- 
ear hyperbolic Klein—Gordon equation (17) in the interval 0 < a < 1, 
a=-1,6 =0,7 = 1,andk = 2. The initial conditions are given by 


Co ie 


and f(z, t) = 6xt(x? —t?)+2°t® in the interval 0 < x < 1. The exact solution 
of this example is 
u(x,t) = 293. (34) 


The Dirichlet boundary conditions are extracted from the exact solution eas- 
ily. 


Table 2: Error norms of Example 2 at different times T using N = 20 grid points, 
dt = 0.001, MQQI-PS and c = 0.381 

MQQI— PS TMQQI 

Lx Ta RMS Loo[37] La[37| RMS[37| 
4.03996 x 10-& ~—-1.73964 x 10-& 1.69559 x 10-© 7.7958 x 10-© = 3.4694 x 10-> 4.8581 x 10-© 
1.76028 x 10-5 7.79986 x 10-® 7.60237 x 10-® =: 1.2307 x 10-4 55.5475 x 10-4 7.7680 x 107 
1.89709 x 10-* 7.44319 x 10-5 7.25472 x 107° 5.3019 x 10-4 2.4618 x 10-3—- 2.4473 x 1074 
5.75087 x 10-* 2.18034 x 10-4 2.12514 x 10-4 1.8601 x 10-3 2.4618 x 10-* 9.7029 x 10-4 
1.5485 x 1073 6.47638 x 10-4 6.31239 x 10-4 3.5192 x 10-8 )=—-1.3585 x 10-2. 1.9022 x 1073 


aoRWN FH) YW 
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0.0015 


o 
So 
co 
° 


Absolute Error 
o 
o 
S 
uo 


Figure 2: Numerical and exact solutions (left)of Example 2 and error of method (right) 
with N = 20, dt = 0.001 at T=1,2,---,5 


Table 3: Error norms of Example 3 at different times T using N = 50 grid points, 
é6t = 0.001, MQQI-PS, and c = 0.181 
T Les Ly RMS Loo [37] L,[37] RMS [37] 
1 2.83572 x 10-* 1.59049 x 10-*_—1.11334x 10-7 4.2672 x 10-® 2.8456 x 10-7: 2.0242 x 10-® 
2 7.81595 x 10-4 3.84643 x 10-4 2.6925 x 10-4 44.9397 x 10-4 2.6266 x 10-3 1.8577 x 10-4 
3 6.78361 x 10-4 1.61686 x 10-* 1.13818 x 1074 —. 2.0880 x 10-3 9.5859 x 107-3 6.7016 x 1074 
: 


1.09944 x 10-4 7.14267 x 10-5 4.99987 x 10-5 6.4162 x 10-8 2.3457 x 10-2: 1.8844 x 107 
3.22279 x 10-4 1.60303 x 10-4 1.12212 x 10-4 1.7941 x 10>. 7.4214 x 10-* 5.2391 x 10-3 


Table 2 presents L., errors, D2 errors, and RMS errors in Example 2. The 
results are obtained by using N = 20 grid points, dt = 0.001, MQQI-PS, and 
c = 0.381. We compare our results with the reported results in [37]. As can 
be seen from the tables, the numerical results are in good agreement with the 
exact solution. Figure 2 compares numerical and exact solutions of Example 
2 with N = 20 grid points and dt = 0.001 at different times T’. 


Example 3. In this example, we consider the one-dimensional nonlinear 
hyperbolic Klein—Gordon equation (17) in « € [—1,1], with a = -1,8 = 
1,y = land k = 3. The initial conditions are given by 


u(z,0) = z* coshz,-1<2 <1, 
u;(,0) = x*sinhz,-1<<2 <1, 


and the right-hand side function is 
f (x,t) = (a? — 2) cosh(a + t) — 4a sinh(x + t) + 2° cosh? (x + t). 
The exact solution is given as 


u(a,t) = 2* cosh(x + t). (35) 


The Dirichlet boundary conditions are extracted from the exact solution eas- 
ily. 

Table 3 presents DL. errors, Lz errors, and RMS errors in Example 3. The 
results are obtained by using N = 50 grid points, dt = 0.001, MQQI-PS, and 
c = 0.181. We compare our results with the reported results in [37]. As can 
be seen from the tables, the numerical results are in good agreement with 
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—t=3 
— t=4 


Percentage Error 


Figure 3: Numerical and exact solutions of Example 3 and error of method with N = 50, 
dt = 0.001 


the exact solution. Figure 3 compares the numerical and exact solutions of 
Example 3 with N = 50 grid points and dt = 0.001 at different times T’. 


Example 4. In this example, we consider the one dimensional nonlinear 
hyperbolic Klein—Gordon equation (17) in x € [0,1], with a = —2.5,8 = 
1,y = 1.5, and k = 3 and f(x,t) =0. The exact solution is given by: 


u(az,t) = Btan(K(x£ + Ct)), (36) 


where B= \/2, K = /7a4en. 


The initial and boundary conditions are extracted from the exact solution 
easily. 

Table 4 presents D., errors, D2 errors, and RMS errors in Example 4 with 
C' = 0.05. The results are obtained by using N = 20 grid points, dt = 0.001, 
MQQI-PS, and c = 2.2. We compare our results with the reported results in 
[37]. As can be seen from the tables, the numerical results are in good agree- 
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Table 4: Error norms of Example 4 with C = 0.05 at different times T using N = 20 
grid points, dt = 0.001, a = —2, 8 = 1, y = 1.5, MQQI-PS method and MQQI shape 
parameter c = 2.2 


Tr 


Loo 


D2 


RMS 


Tx [37] 


L2[37] 


RMS [37] 


1 


2 
3 
4 


8.10099 x 10-7 
2.34139 x 10-° 
2.77812 x 10-6 
2.63626 x 10-6 


4.07215 x 10-7 
1.21136 x 10-® 
1.5571 x 10-® 

1.21796 x 10-® 


3.96904 x 10-7 
1.18068 x 10~® 
1.51767 x 10-6 
1.18712 x 10-® 


2.1781 x 10-7 
3.0648 x 10-7 
3.7008 x 10-7 
3.4230 x 1077 


1.2850 x 10-® 
1.4099 x 10- 
1.3144 x 10-6 
2.0012 x 10~§ 


1.2787 x 10-7 
1.3980 x 1077 
1.7004 x 10-7 
1.2002 x 10-7 


0.5F 4 


0.45 1 


0.3} | 


0.2} | 


G1; 1 


0.0} zl 


2.5% 10-6 
2.x 10-6 


1.5x 10-8 
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Figure 4: Numerical and exact solutions of Example 4 and error of method with C=0.05, 
N = 20 and dt = 0.001 at different times T 


ment with the exact solution. Figure 4 compares the numerical and exact 
solutions of Example 4 with N = 20 grid points and dt = 0.001 at different 
times T’. 

Table 5 presents LD. errors, [2 errors, and RMS errors in Example 4 with 
C= 0.5. The results are obtained by using N = 50 grid points, dt = 0.001, 
MQQI-PS method and c = 0.97. We compare our results with the reported 
results in [37]. As can be seen from the tables, the numerical results are in 
good agreement with the exact solution. Figure 5 compares the numerical 
and exact solutions of Example 4 with C = 0.5, N = 50 grid points and 
dt = 0.001 at different times T. 
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Table 5: Error norms of Example 4 with C = 0.5 at different times T using N = 50 
grid points, dt = 0.001, a = —2, 8 = 1, y = 1.5, MQQI-PS method and MQQI shape 
parameter c = 0.97 
T Le Ty RMS Lo[37] L [37] RMS [37] 
1 1.34058 x 10-® 3.6230 x 10-” 7.95305 x 10-” 5.2134 x 10-© 4.0608 x 10-5 4.0392 x 10° 
2 5.28169 10-® 1.462 x 10-6 1.92792 x 10-© 2.1805 x 107 1.5588 x 10-4 1.5505 x 107° 
3 2.54257x 107° 5.6233 x 10-& 4.46219x10-® 9.0113 x 10-5 6.3231 x 10-4 6.5903 x 10-5 
4 1.32595 x 1072 2.53105 x 10-* 9.82629 x 1075 8.2372 x 10-4 5.3386 x 107% 5.3219 x 107-4 


0.000025 
0.00002: 
0.000015: 


0.00001 — tx2 


Absolute Error 


5.x 10°) 


—t=4 


Absolute Error 
o 
S 
S 
a 


Figure 5: Numerical and exact solutions of Example 4 and error of method with C=0.5, 
N = 50 and dt = 0.001 at different times T 


Example 5. Let us consider the following two-dimensional nonlinear hy- 
perbolic Klein—Gordon equation (17) in (z,y) € [0,7] x [0,7], with a = 
—1,8 =—-2,y =1, and k = 2. The exact solution of this example is 


u(x, y,t) = sin(x) sin(y) cosh(t), (37) 
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Table 6: Error norms of Example 5 at different times T using N = 10 x 10 uniform grid 
points, dt = 0.01, a 1,2 2,7 = 1, k = 2, MQQI-PS method and MQQI shape 
parameter c = 4.15 


MQQI-PS  (N=10x 10) MK—MLPG (N =33 x 33) 
T Leo Ly RMS L[41] L,[41] RMS [41] 
1 2.32242 x 10-5 1.4976 x 10-® ~— 8.40697 x 10-7 1.04 x 10-3 1.28 x 10? 3.87 x 10-4 
2 3.88411 x10-® 2.65018 x10-® 1.48771 x 10-® 1.23 x 10-? 1.86 x 107! 5.64 x 107% 
3 4.19091 x 10-® 3.26018 x10-® 1.83014x10-® 4.35 x 107? 7.52 x 1071 2.28 x 10-? 
4 6.10456x10-® 4.66624 x 10-® 2.61945 x 10-8 1.20 x 1074 2.07 x 10° 6.27 x 107? 
5 9.83839 x 10-® 6.92622 10-® 3.88811x 10-6 3.25 x 107! 5.77 x 10° 1.75 x 1071 


the initial and Dirichlet boundary conditions are extracted from the exact 
solution easily. 
Table 6 presents D,., errors, £2 errors, and RMS errors in Example 5. The 


Numerical solution(T=5) Exact solution(T=5 


Figure 6: Numerical and exact solution of Example 5 with N = 10 x 10, dt = 0.01, 
a 1,8 2,7 =1, k = 2, and MQQI shape parameter c = 4.15 at time t=5 


results are obtained by using N = 10 x 10 grid points, dt = 0.01, MQQI-PS 
method, and MQQI shape parameter c = 4.15. We compare our results with 
the reported results in [41]. As can be seen from the tables, the numerical 
results are in good agreement with the exact solution. Figure 6 compares the 
numerical and exact solutions of Example 5 with N = 10 x 10 grid points 
and dt = 0.01 at ¢ = 5.The contour plot of the absolute error norm for the 
MQQI-PS method in Example 5 is portrayed in Figure 7. The MQQI shape 
parameter c versus L. is shown in Figure 8 when t = 5, dé = 0.01, and 
N = 10 x 10, corresponding to Example 5. 


Table 7 presents comparisons of D,, errors, D2 errors, and RMS errors for 
the MQ-PS method and MQQI-PS method in Example 5. The results are 
obtained by using N = 25 x 25 grid points, dt = 0.01, MQ-RBF-PS method 
with shape parameter c = 0.37 and MQQI-PS method with shape parameter 
c= 4.6. As can be seen from the tables, the numerical results from MQ-PS 
and MQQI-PS are in good agreement with the exact solution. The contour 
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Figure 7: The contour plot of the absolute error norm for the MQQI-PS method in 
Example 5 with N = 10 x 10, dt = 0.01 and MQQI shape parameter c = 4.15 at time t=5 
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Figure 8: The plot of Leo error norm against MQQI-Ps(c) in Example 5 with N = 
10 x 10, dt = 0.01 and MQQI shape parameter c = 4.15 at time t=5 


plot of the absolute error norm for the MQ-PS method and MQQI-PS method 
in Example 5 are compared in Figure 9 by using N = 25 x 25 grid points. 


5.2 Sine-Gordon equation 


Consider the hyperbolic Klein—Gordon equation () with nonlinear term 
(wu) = sin(u). In this section, we consider Sine-Gordon equations of the 


form 
ute(z,t) + aAu = —sin(u). (38) 


With the initial conditions 
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Table 7: Compare error of MQ-PS method and MQQI-PS method for Example 5 with 
a=-—1,$6=-—2,7=1, and k = 2 at different times T using N = 25 x 25 grid points and 


dt = 0.01 
MQ-PS (c = 0.37) MQQI-PS  (c=46) 
T Diss Te RMS Les Ty RMS 
1 9.99276 x 10-® 4.48038 x 10-& 2.52576 x 10-® 1.15822 x 10-® 5.21147 x 10-® ~—-2.93791 x 10-* 
2 6.77584 10-5 2.52980 x 10-5 1.42615 x 10-® 7.47466 x 10-5 3.1979 x 10-® ~—-:1.80278 x 10-® 
3 2.78783 x 10-4 7.76853 x 10-5 4.37942 x 10-5 2.54482 x 10-4 1.0254 x 10-4 5.78057 x 1075 
4 5.77975 x 107-4 1.50924 10-4 8.50817 x 10-5 5.74854 x 10-4 1.47003 x 10-4 8.28710 x 10-® 
5 5.72217 x 1074 1.24466 x 10-4 7.01665 x 10-5 7.96027 x 10-4 1.64803 x 10-4 9.29055 x 1075 


0.0007 


0.0001 


Figure 9: The contour plot of the absolute error norm for the MQ-PS method with 
MQ shape parameter c = 0.37 (up) and MQQI-PS method with MQQI shape parameter 
c = 4.6 (down) in Example 5 with N = 25 x 25, dt = 0.01 at time t=5 


u(x,0) = u(r) «x E€NCR4(d=1,2), 


(39) 
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Table 8: Error norms of Example 6 with C = 0.5 at different times T using N = 121 grid 
points, dt = 0.001, and MQQI shape parameter c = 1.275 
T Loo I RMS L[43] Lx[9] 
9 3.50847 x 10-4 5.1875 x 10-4 6.66931 x 10-° 6.296 x 10-4 4.518 x 107! 
36 2.4275 x 10-4 8.69372 x 10-4 1.11771 x 10-4 8.109 x 10-4 2.293 x 10° 
1083.60803 x 10-4 1.5558 x 1073 2.00022 x 10-4 2.568 x 10- ~—5.120 x 10° 


uz(z,0) =up(x) rE, (40) 


and Dirichlet boundary conditions 
u(a,t) = g(t), t€(0,T], x ean. (41) 


We can transform problem (38)-—(41) into a system of coupled first order 
equations in time given by 


eS = v(z,t), (42) 


vz(x, t) = —aAu — sin(u), 
with the initial conditions 
u(x,0)=up(t) re QcRy’, (43) 


v(@,0) =v(z) rEQ, (44) 


and boundary conditions 
u(xz,t)=9(t), t€(0,T], re an. (45) 


The RBF-PS method described at the previous section is very efficient and 
can be applied to approximate solution of the system (42)-(45). It is enough 
to replace the vector U* in (30) by the vector sin(U”). Numerical solutions 
obtained by this method are compared with the exact solutions and other 
existing methods. 


Example 6. Consider problem (38)—(41) with interval « € [—2,58] and 
a=-—l. 
The exact solution is given by: 


u(x,t) = 4tan71(exp(y(x — Ct) + B)), (46) 


_ 1 
where 7 = Jia 
The initial and boundary conditions are extracted from the exact solution 


easily. 

Table 8 presents L,, errors, [2 errors, and RMS errors in Example 6 with 
C' = 0.5. The results are obtained by using N = 121 grid points, dt = 0.001, 
and the shape parameter c = 1.275. We compare our results with the reported 
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Figure 10: | Numerical and exact solutions of Example 6 (up) and error of method 
(down) with C = 0.5, N = 121 and dt = 0.001 


results in [43, 9]. As can be seen from the tables, the numerical results are in 
good agreement with the exact solution. Figure 10 represents numerical and 
exact solutions of Example 6 and error of method with N = 121, dt = 0.001 
at different times T’. 


Example 7. In this example, we present numerical results of MQ-RBF- 
PS method for the two-dimensional sine-Gordon equation. Consider problem 
(38)-(41) with domain (x,y) € [—7, 7] x [—7,7] and a = -1. 

The exact solution is given by [24] 


u(z,y,t) = 4tan7(exp((x + y — t))). (47) 


The initial and Dirichlet boundary conditions are extracted from the exact 
solution easily. 

Table 9 presents Lo, errors, D2 errors, and RMS errors in Example 7. 
The results are obtained by using N = 25 x 25 grid points, dt = 0.001, 
and the shape parameter c = 0.31. As can be seen from the tables, the 
numerical results are in good agreement with the exact solution. The graphs 
of the numerical and analytical solutions at t = 5 are shown in Figure 11. 
For Example 7, the contour plot of the absolute error norm of the MQ-PS 
method is portrayed in Figure 12. The MQ shape parameter c versus Loo 
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Table 9: Error norms of Example 7 at different times T using N = 25 x 25 grid points, 
dt = 0.001 and MQ shape parameter c = 0.31 


T 


Lin 


L2 


RMS 


oR WN 


1.86749 x 107 = 8.23484 x 107-4 
5.07197 x 10-3 —- 1.92593 x 10-8 
5.48379 x 10-3 2.67973 x 107-3 
4.88027 x 1072 =. 2.97766 x 1073 
5.78486 x 10-3 —- 3.21738 x 107 


2.19909 x 1074 
5.14315 x 10-4 
7.15615 x 107-4 
7.95177 x 107-4 
8.59192 x 10-4 


olution (T=5 


Figure 11: Numerical and exact solution of Example 7 with N = 25 x 25, dt = 0.001 
and MQ shape parameter c = 0.31 at time t=5 


Absolute errorMO-PS method, T=5) 


Figure 12: The contour plot of the absolute error norm of the MQ-PS method with 
N = 25 x 25, dt = 0.001 and MQ shape parameter c = 0.31 at time t=5 for Example 7 


corresponding to Example 7 is shown in Figure 13 when t = 5, dt = 0.001, 
and N = 25 x 25. 
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MQ shape parameter c 


L,, Error 


0.05} | 


Figure 13: The plot of Loo error norm against MQ-Ps(c) with N = 25 x 25, dt = 0.001 
at time t=5 for Example 7 


6 Conclusion 


In this paper, a new radial basis function-pseudospectral (RBF-PS) method 
was proposed to solve the second order Klein—Gordon equation with quadratic 
and cubic nonlinearity and sine-Gordon equation. In addition, a new tech- 
nique was introduced to force approximations to satisfy the boundary condi- 
tions exactly each time. The proposed numerical method, which is truly 
meshless, was based on the MQ radial basis function, and MQ quasi- 
interpolation was used to compute RBF differentiation matrices to discretize 
the spatial derivative operators. The fourth-order Runge-Kutta time step- 
ping method was used to deal with the temporal part of the problem. Some 
illustrative examples have been given to show the accuracy and applicability 
of the presented method. 
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